#include "Init.h"
#include "Basic_Opera.h"

/**
* @Name: Init 
  1. Initialize u, v and omega in C
  2. Initialize rho in R
*/
void Init()
{
	double X, Y, tmp;	

	for(k = 0; k < M*NUM; k++){
		C_omega[k][0] = 0.0;
		C_omega[k][1] = 0.0;
		C_u[k][0]     = 0.0;
		C_u[k][1]     = 0.0;
		C_v[k][0]     = 0.0;
		C_v[k][1]     = 0.0;
	}

	double *rho1, *rho2;
	rho1 = (double *)malloc(sizeof(double)*(M*NUM));
	rho2 = (double *)malloc(sizeof(double)*(M*NUM));

	for(j = 0; j < M; j++){
		X = (j + M*myrank)*DH;
		for(i = 0; i < NUM; i++){
			Y   = i*DH;
			tmp = pow(PI,2) - pow(Y-PI,2) - pow(X,2);
			if(tmp > 0){
				rho1[index] = exp(1-PI*PI/tmp);	
			}
			else{
				rho1[index] = 0;
			}
		}
	}

	for(j = 0; j < M; j++){
		X = (j + M*myrank)*DH;
		for(i = 0; i < NUM; i++){
			tmp = pow(1.95*PI,2) - pow(X-2*PI,2);
			if(tmp > 0){
				rho2[index] = exp(1-pow(1.95*PI,2)/tmp);	
			}
			else{
				rho2[index] = 0;
			}
		}
	}
	
	for(j = 0; j < M; j++){
		for(i = 0; i < NUM; i++){
			R_rho[index][0]	= 50*rho1[index]*rho2[index]*(1-rho1[index]);
			R_rho[index][1] = 0;
		}	
	}

	FFT_R2C_For_Input(R_rho,C_rho);               /// transpose rho from R to C
	
	free(rho1);
	free(rho2);
}




